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Abstract 

In this paper we compute the Fourier spectrum of the 
Fractal Interpolation Functions FIFs as introduced 
by Michael Barnsley. We show that there is an an- 
alytical way to compute them. In this paper we at- 
tempt to solve the inverse problem of FIF by using 
the spectrum 

1 Iterated Function Systems 

The affine transform performs translation stretching 
and rotation on a given set. In the special case of 
two dimensions the afhne transform on a set S in 
2-D space is described by the equation: 
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where (x, y) G S. The effects of an affine 
transform on a set are depicted in fig. [1] The 
union of N affine transformations is called the 
Hutchinson operator: W = U^=i^n- For a 
specified metric the distance /i(.A, B) between two 
sets B can be defined. Under certain condi- 
tions [2] the Hutchinson operator is contractive, 
h{WiA),WiB)) < sh{A,B), s < 1. Successive it- 
erations with Hutchinson operator on a random set 
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results in a sequence of a sets that converges in the 
attractor of the operator A, which satisfies the con- 
dition A — W{A) — {Jn=i Wn{A). Any system that 
uses the Hutchinson operator in order to generate it- 
eratively the attractor A is called Iterated Function 
System (IFS). 



Original /V Translation & Scaling 

Translation & Rotation 

Figure 1: An affine transform can translate rotate 
and flip a 2-dimensional shape. 



1.1 Fractal Interpolation Functions 

Fractal Interpolation Functions (FIF) is a special case 
of the 2-dimensional IFS and maintain all their char- 
acteristics. FIF attractors are continuous functions 
that can be used to model continuous signals. FIF 
interpolate a given set of -f 1 points (x„, ?/„), n = 
0,1... Af 

Xo < Xi < X2 < X3 < . . . < XN 
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The FIF that interpolates the above set is comprised 
of N affine maps: 
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The necessary condition is that the Interpolation 
Function passes from the N + 1 initial points, 



(1) 



Xq 






and Wn 


XN 










Vn-l 


VN 




Vn 



n = 1, 



We call order of the FIF, the number N of the 
affine maps. The conditions provide 4 equations for 
5 parameters , so (i„, the vertical scaling factor is 
chosen to be the free parameter. If we solve the above 
equations for a„, c„, e^, /„ in terms of d„, we find : 



satisfies the conditions of ([T]), 
has a unique fixed point, the function F £ C, 
{TF){x) ^ F{x), yxe[xo,XN]. 

The following restrictions guarantee the contractivity 
of T operator: 

• a„ < 1 

• \d„\ < 1. 

Barnsley's operator is very useful because the at- 
tractor of an FIF can be generated with the iter- 
ative application of T on an initial signal F(0) : 

F^ lim T°"-'{Fa),Fa eC 
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Let the real numbers an, Cn, Cn, fn be defined by 
(l2][5]). Barnsley [2' introduced the operator T for the 
class C of continuous functions, T : C — > C by 

iTF)ix) = Cn£~\x) + dnF{ln\x)) + /„ 

Vx e [x„_i,a;„] ,n= l,2...,iV, 

and 

In ■■ [xq,Xn] [Xn~l,Xn] 

the invertible transformation 

^nix^ — dnX 6n- 

The above operator: 

• is a contraction mapping according to Hausdorff 
metrics. 




Figure 2: Formation of the FIF attractor after suc- 
cessive iterations 



1.2 Discrete FIF 

One of the most interesting properties of the FIF is 
that its attractor is independent of the initial signal 
(initiator). If the initiator is a continuous signal, for 
example the linear interpolation between the given 
interpolation points, all the instances throughout all 
the steps of the iterations will be continuous signals, 
fig- m On the other hand if the initiator is a sin- 
gle point the attractor formed after infinite iterations 
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will be a continuous signal, but signals instances of 
the FIF throughout the iterations will be discrete sig- 
nals, fig. O Although during all the iterations the 
instances are discrete signals or strictly mathemati- 
cally speaking finite countable sets, the attractor is 
a continuous signal, an infinite and uncountable set. 
Instances of a the formation of an FIF from a discrete 
initiator are shown in fig. [3) 

2| ^ ^ ^ 




0.26 1.57 2.87 



Figure 3: Generation of a discrete FIF. 

It is essential to show that a good choice of the 
initiator is the + 1 given interpolation points 
{xn,yn),n' = 0,1... iV. After the first iteration in 
each of the A^ subintervals between the points, N — 2 
new points are generated fig. [S] Let these points be 
{xs,ys), s = 1, . . . ,{N - 1)N. Notice that aU these 
A^ + 1 + (A^ - 2) A^ ^ N'^-N+l points belong to the 
attractor of FIF. This wouldn't be true if the initiator 
was a different set apart from the given interpolation 
points. By repeating this procedure after m itera- 
tions we get an M(m)-point discrete sequence that 
is a sampling of the FIF's attractor, with sampling 
period Tg = jj. If the initiator included any other 
irrelevant point, that would not be mapped to an at- 
tractor's point after a finite number of iterations. It 
can be proved that if the initiator is not the set of 
the initial interpolation points, then the error of the 



formed sequence after the mth iteration from the at- 
tractor decreases and goes to zero as m ^ cx). The 
number of the attractor's samples is : 

A^ + 1 

A^ + 1 + (A^ - 1)A^ = + 1 
N +1 + {N -1)N + {N -1)N^ = + 1 

M = M{m) = A^" + 1. (6) 

The discrete signal formed after the mth iteration 
over the discrete initiator, is called discrete FIF : 

J[n]=F{xn), n = 0,l,...,M (7) 

Before expanding Barnsley operator T introduced 
above for the discrete FIF it is necessary to make 
clear that according to the strict definition of Frac- 
tals the discrete FIFs are not fractal sets because they 
are finite. Discrete FIFs must be considered as ap- 
proximations of the continuous. 

Since the Barnsley operator assumes infinite res- 
olution, it doesn't apply in discrete signals. So for 
discrete FIF the following modified operator is used: 

N 

r/N = ^(4T^/[^^^] +c.(^) + /„) (8) 

k=l ^ 

(u(n - Cfc-i) - u{n - Cfc)). 

The Symbol ]^ f[n\ defines that within two suc- 
cessive samples of / [n] , A^ — 1 zeros have been inter- 
polated, fig. m The term /„ denotes the parameter 
of FIF as defined in dS]) and should not be confused 
with the discrete signal f[n]. 

Similarly to the continuous discrete FIF of 

m iterations can be constructed with the following 
procedure: 

/m+lW = T/mW- (9) 

It is obvious that when m ^ oo the discrete FIF 
becomes continuous. Notice that T depends on the 
number of iterations m. More specifically the pa- 
rameters of FIF as defined in ([SJS]) depend on the 
Xn ■ The Xn change value because of the upsampling. 
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• Fo = Fn = 0. 

• and the points are evenly spaced. 
Then (p]|5)) becomes, 
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Figure 4: Barnsley operator for discrete signals. The 
initiator is upsampled and interpolated according to 



The generation of an FIF can be represented in terms 
of a linear system as shown in fig. [5] 
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The above equations show that the FIF parameters 
are decoupled from each other. This is very impor- 
lOtant because they can be estimated independently. 
Moreover it is clear that if a„ parameter is estimated 
then all parameters can be found directly, except for 

Any FIF can be transformed to an equivalent FIF 
that satisfies the first two conditions without loosing 
its fundamental properties. More specifically fig. [6] 
shows how a given FIF can be transformed so as to 
satisfy the conditions for the first and the last point. 
It is convenient to use an auxiliary affine map Waux 
that will rotate, scale and translate the given one. 
The Waux transform is invertible and does not affect 
the intrinsic parameters an,dn. 



Figure 5: Block diagram of the FIF Genaration. 



2 Computation of 
Fourier Spectrum 



FIF's 



TIME (sec) 



In order to simplify (l2][5]), it is very convenient to 
adopt the following assumptions. 

• = 0, xn — 1, 



Figure 6: An example of FIF out of range [0,1]. 
Applying an affine transform we can tie it at points 
(0,0) and (1,0). 
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Let A — {{x,y) : F{x) — y}, be the original at- 
tractor and 



and 



u{x) 



X <0 

1 x>0 



(14) 



where Wn, 



G{x) is the piecewise linear function between the in- 
terpolation points : 
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,n = 1, . . . ,N ^(^) = X! + fn)iu{x - Cn-i) - u{x - en)) ■ 



The new transformed attractor is A' = 
{{x' ,y'), F' (x) = y'}. The new points are con- 
nected with the initial 
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That means A' — WauxA. These new points belong to 
a new FIF. Setting in ([H]) A = w^^^A' and applying 

the map Waux, 



M WauxWnwJ^iA') = A', 
^-^n— 1 



(15) 



with 





1/ daux 



By setting faux = Saux = 0, daux = 1, the afhne maps 
of the new FIF are : 



a„ 

(Cn dnCaux ~t- C-aux^n) I ^aux dn 
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Notice that the new FIF has the same order and the 
same d„ parameters. 

2.1 Spectrum of Continuous FIF 

The application of the above simplifications to the 
Barnsley operator T results in the following equation 
for the FIF: 




FREQUENCY (Hz) 

Figure 7: An example of Q(ri). 
Applying continuous fourier transform [5]: 

/ + 00 
F{x)e-'^''dx 
-OO 



F{x) < — > T{n) 

T{n) = g{n)+ aT{an) ^ d„e-*""" , (17) 



N 



n=l 



(18) 



We define the function 



F{^)-En=iidnF{^-^)+ Qin)^ay^dne 
+c„(2^) + /„)(w(a; - e„_i) - it(x - e„)) 



N 



-in 



N-l 

^=a^d„+ie-'^ (19) 

n=0 
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Notice that the above function is the discrete 
time fourier transform of the discrete sequence 
{di,d2, . . . d]\[}, so we deduce that Q{il — 2kTrN) = 
Q{2{k + l)TrN -Q), fc = l,2,... (fig.[71). 




° TIME (aec^'^ " ^ ' ° ™" FREQUENCY (Hj^ 



Figure 8: Left column: Fractal interpolation between 
points (0, 0), (0.25, 1), (0.5, 1.4), (0.75, -0.5), (1,0) af- 
ter 1,3,5 iterations, right column: Corresponding 
spcctrums. 

The Fourier spectrum satisfies the following equa- 
tion: 

= Q{n)T{an) + g{n) (20) 

Through Barsley Operator in frequency domain 
the Fourier spectrum can be computed iteratively, 
fig.l 

^,n+iin) = Q{n)J'„^ian) + Gin) (21) 

^o(^^) = 0. 

After infinite iterations: 

00 2 — 1 

J^{n) = J2 Oi^a') Y[ Q{na') (22) 

1=0 j=0 



2.2 Spectrum of Discrete FIF 

From the ([H]) it is obvious that the FIF signals are 
not band-limited. As a result the spectrum of a dis- 
crete FIF is aliased. Assume Ts = jj has been cho- 
sen as the sampling period, then the Discrete Time 
Fourier Transform (DTFT) 5^ is : 

M 

/>)-^/[n]e~^"",c^ = f]r,. 

As in continuous case for the computation of the 
DTFT the T operator is used in frequency domain: 

N 

/m+i(t^) = frn{ai^) ^ dkae^'"^"" + .9m(a;) (23) 
p=i 

Posing the analogy between the continuous and the 
discrete case we define the q function : 

N 

q^{uj) = a^dpe'"^''''. 
p=i 

N 

q{uj) = Y,dpe-'^''', p = 0,l,...,iV. 
and Bp = 2^(Af — 1), so 

- E^-o' d,-,ie(-^)(^"-^-i:::r ^ = o, i, . . . , a^. 

(24) 

/m+l(w) = qra{l^)fm{at^) + gm{io) (25) 

After m iterations 

m i — 1 

Jni{uj) = ^g™_j(tja*) gi_j+i(tja-') (26) 

i=0 i=0 

The Discrete Fourier Transform can be com- 
puted after keeping the frequency values between 
and 27r and sampling the spectrum at cj = 
0,f,2|f,...,(Af-l)f ,fig.[Sl 

In the continuous case Q{Vl) has period 2nN . In 
the discrete case where uj G [0, 2tt] the q{uj) function 
has period 27r . . 
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of linear equations. 



FREQUENCY (H!) 



Figure 9: 5-point FIF attractor with dn = 
-0.74,0.8,-0.77,0.85,0.88. and the spectrum after 
the 5th iteration. 



3 FIF parameter estimation us- 
ing spectrum 



FIF modeUing of signals has been proposed by Mazel 
[3], using information of signal in time domain. In 
this section the spectrum of signal is used to estimate 
its parameters, provided the FIF order is known. 

Given the signal f[n],n = 0, 1, . . . , M , let be 
FIF's estimated order. We expect its spectrum to 
satisfy the following equation: 



aujk 
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{MN - 1) 
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All values that zero denominator are excluded. No- 
tice that fm+i[n] and fm[n\ signals must have the 
same length. Considering also that /m+i[0] — 
/„[0] and fm+i[M{m + 1)] = /™[M(to)], where 
M{m),M{m + 1) their lengths. In order to equalize 
their lengths it is necessary to interpolate with zeros 
fm [n] ■ Because of the term /,„ [awk\ the computation 
of FFT must be on ||^, A: = 0, ... , NM - 1. This 
is done by padding the signal with MN — M zeros. 
The dn are determined by the solution of the system 
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(c) 
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Figure 10: (a) FIF without noise. In (b),(c) FIF 
contaminated with noise SNR — 20 and SNR — 10. 
In (d),(e) reconstructed FIF. 



4 Results 

The above algorithm was tested in FIF signals con- 
taminated with white noise and the parameters ex- 
tracted were very close to the real, fig. [TUl The main 
advantage of the method is the use of FFT. It is well 
known that FFT is quite simple and easily imple- 
mented. The first disadvantage of the method is that 
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it cannot find the FIF order. In the above experi- 
ments we tried for = 2, 3, . . . until the estimated q 
function satisfied the periodicity and symmetry con- 
ditions mentioned earher. The second and most im- 
portant problem is that it is very sensitive in the win- 
dow effect of the FFT. It is known that although a 
part of FIF signal is self affine it is not an FIF. Trying 
to model it as an FIF results in wrong estimations. 
In the example of fig. [Til it is evident that the period 
of the q function has been expanded. Having chosen 
the order N, the period of q is known 2Tr . In 
the right plot of fig. [TT] it is evident that the period 
of q is much higher than the expected. But although 
the spectrum is not the best method for solving the 
inverse problem it can be used as powerful analysis 
tool. As shown above it can reveal the FIF nature 
of a signal and it can also help in the prediction of a 
missing part of a time series, given that it belongs to 
class of FIF. 




Figure 11: Left : q{uj). Right : of the corrupted 
FIF. 
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